AMR: Resistencia a Antimicrobianos

Librerías empleadas

library(readr)
library(dplyr)

## 
## Adjuntando el paquete: 'dplyr'

## The following objects are masked from 'package:stats':
## 
##     filter, lag

## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union

library(ggplot2)
library(plotly)

## 
## Adjuntando el paquete: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

library(sf)

## Warning: package 'sf' was built under R version 4.4.2

## Linking to GEOS 3.12.2, GDAL 3.9.3, PROJ 9.4.1; sf_use_s2() is TRUE

library(pxR)

## Cargando paquete requerido: stringr

## Cargando paquete requerido: reshape2

## Cargando paquete requerido: RJSONIO

## Cargando paquete requerido: plyr

## ------------------------------------------------------------------------------

## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)

## ------------------------------------------------------------------------------

## 
## Adjuntando el paquete: 'plyr'

## The following objects are masked from 'package:plotly':
## 
##     arrange, mutate, rename, summarise

## The following objects are masked from 'package:dplyr':
## 
##     arrange, count, desc, failwith, id, mutate, rename, summarise,
##     summarize

library(RColorBrewer)
library(httr)

## 
## Adjuntando el paquete: 'httr'

## The following object is masked from 'package:plotly':
## 
##     config

library(tidyverse)

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ forcats   1.0.0     ✔ tibble    3.2.1
## ✔ lubridate 1.9.3     ✔ tidyr     1.3.1
## ✔ purrr     1.0.2

## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ plyr::arrange()   masks plotly::arrange(), dplyr::arrange()
## ✖ purrr::compact()  masks plyr::compact()
## ✖ httr::config()    masks plotly::config()
## ✖ plyr::count()     masks dplyr::count()
## ✖ plyr::desc()      masks dplyr::desc()
## ✖ plyr::failwith()  masks dplyr::failwith()
## ✖ plotly::filter()  masks dplyr::filter(), stats::filter()
## ✖ plyr::id()        masks dplyr::id()
## ✖ dplyr::lag()      masks stats::lag()
## ✖ plyr::mutate()    masks plotly::mutate(), dplyr::mutate()
## ✖ plyr::rename()    masks plotly::rename(), dplyr::rename()
## ✖ plyr::summarise() masks plotly::summarise(), dplyr::summarise()
## ✖ plyr::summarize() masks dplyr::summarize()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors

library(jsonlite)

## 
## Adjuntando el paquete: 'jsonlite'
## 
## The following object is masked from 'package:purrr':
## 
##     flatten
## 
## The following objects are masked from 'package:RJSONIO':
## 
##     fromJSON, toJSON

library(readxl)
library(rjson)

## 
## Adjuntando el paquete: 'rjson'
## 
## The following objects are masked from 'package:jsonlite':
## 
##     fromJSON, toJSON
## 
## The following objects are masked from 'package:RJSONIO':
## 
##     fromJSON, toJSON

library(tidyr)

Cargar environment

load("OUTPUT/enviroment.RData")

Pregunta 1: Consumo de antibióticos en España

Determinar el consumo de antibióticos por comunidades autónomas, diferenciando entre recetado y no recetado, para identificar tendencias de uso.

tipo_ccaa_consumo_o_no <- read_delim("INPUT/DATA/datos_ccaa/tipo_ccaa_consumo_o_no.csv", 
                                     delim = ";", escape_double = FALSE, trim_ws = TRUE)

## Rows: 3420 Columns: 6
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ";"
## chr (6): Sexo, Total Nacional, Comunidades y Ciudades Autónomas, Tipo de med...
## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

antibioticos <- tipo_ccaa_consumo_o_no %>% 
  filter(`Tipo de medicamento` == "Antibióticos")

antibioticos <- tipo_ccaa_consumo_o_no[tipo_ccaa_consumo_o_no$`Tipo de medicamento` == "Antibióticos" & tipo_ccaa_consumo_o_no$`Sexo` == "Ambos sexos",]

# Esta tabla enseña que el 3,54% de la población española en la última encuesta reconoce 
# haber consumido antibióticos en las últimas 2 semanas (En el año 2021)


# Comunidad que más consume

antibioticos_si <- antibioticos[antibioticos$`Sí o no` == "Sí",]

antibioticos_si$Total <- gsub(",", ".", antibioticos_si$Total)

antibioticos_si$Total <- gsub(" ", "", antibioticos_si$Total)

antibioticos_si$Total <- as.numeric(antibioticos_si$Total)

# Revisa si hay NA después de la conversión
sum(is.na(antibioticos_si$Total))  # Muestra el número de NAs

## [1] 0

consumo_comunidades <- antibioticos_si %>%
  group_by(`Comunidades y Ciudades Autónomas`) %>%
  arrange(desc(Total))


c_c_final <- consumo_comunidades %>%
  mutate(`Comunidades y Ciudades Autónomas` = ifelse(is.na(`Comunidades y Ciudades Autónomas`), "Total País", `Comunidades y Ciudades Autónomas`))

Con estos datos podemos observar de manera gráfica la tendencia de consumo de antibióticos por comunidad autónoma en este país.

ggplot(c_c_final, aes(x = reorder(`Comunidades y Ciudades Autónomas`, -Total), y = Total)) +
  geom_bar(stat = "identity", fill = "steelblue") +
  labs(x = "Comunidades Autónomas", y = "Total Consumo", title = "Consumo de antibióticos por Comunidad Autónoma") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) 

A continuación estudiaremos si esa cantidad de antibiótico se ha consumido de forma regulada, o por el contrario, no está recetado, favoreciendo así el autoconsumo de los mismos, lo cual se está considerando como una de las posibles causas de la resistencia a los antibióticos

# Leer el archivo .px
archivo_px <- read.px("INPUT/DATA/datos_ccaa/tipo_ccaa_recetado_o_no.px")
df_px <- as.data.frame(archivo_px)



# quedarme solo con los antibióticos
antibiotic <- df_px[df_px[["Tipo.de.medicamento"]] == "Antibióticos" & df_px$`Sexo` == "Ambos sexos", ]


# qué comunidad autónoma se automedica más
antibioticos_sin_receta <- antibiotic[antibiotic$`Recetado` == "No recetado",]

antibioticos_sin_receta <- antibioticos_sin_receta %>%
  filter(value != 0)


consumo_comunidades <- antibioticos_sin_receta %>%
  group_by(`Comunidad.autónoma`) %>%
  arrange(desc(value))

ggplot(consumo_comunidades, aes(x = "", y = value, fill = `Comunidad.autónoma`)) +
  geom_bar(stat = "identity", width = 1) +
  coord_polar(theta = "y") +
  labs(title = "Consumo de Antibióticos por Comunidad Autónoma",
       fill = "Comunidad Autónoma") +
  theme_void()  # Elimina el fondo y los ejes

Pregunta 2: Relación entre PIB y resistencia a antibióticos en la UE

Analizar si existe una relación significativa entre el Producto Interno Bruto de los países de la Unión Europea y los niveles de resistencia a antibióticos.

pibPP <- read_table("INPUT/DATA/datos_pib.tsv")

## Warning: Missing column names filled in: 'X14' [14]

## 
## ── Column specification ────────────────────────────────────────────────────────
## cols(
##   `freq,na_item,ppp_cat,geo\TIME_PERIOD` = col_character(),
##   `2012` = col_double(),
##   `2013` = col_double(),
##   `2014` = col_double(),
##   `2015` = col_double(),
##   `2016` = col_double(),
##   `2017` = col_double(),
##   `2018` = col_double(),
##   `2019` = col_double(),
##   `2020` = col_double(),
##   `2021` = col_double(),
##   `2022` = col_double(),
##   `2023` = col_double(),
##   X14 = col_logical()
## )

# Cambio de nombre columna
colnames(pibPP)[1] <- "pais"

# Nos quedamos con las últimas letras
pibPP$pais <- substr(pibPP$pais, nchar(pibPP$pais) - 1, nchar(pibPP$pais))

lista_pais <- list("BE", "BG", "CZ", "DK", "DE", "EE", "IE", "EL", "ES", "FR", "HR", "IT", "CY", "LV", "LT", "LU", "HU", "MT", "NL",
               "AT", "PL", "PT", "RO", "SI", "SK", "FI", "SE")

# nos quedamos solo con los países de la UE


pib <- pibPP %>% filter(pais %in% unlist(lista_pais))


# quitar la columna nula

pib <- pib[, colSums(is.na(pib)) < nrow(pib)]

pib_2022 <- pib %>% select(pais, `2022`)


# ggplot del pib en el 2022 (primero lo ponemos en descendente)

pib_2022_desc <- pib_2022 %>% arrange(desc(`2022`))


# sustituir las etiquetas de los países

pib_2022_desc <- pib_2022_desc %>%
  mutate(pais = case_when(
    pais == "SK" ~ "Eslovaquia",
    pais == "SI" ~ "Slovenia",
    pais == "EE" ~ "Estonia",
    pais == "MT" ~ "Malta",
    pais == "LV" ~ "Latvia",
    pais == "HR" ~ "Croatia",
    pais == "EL" ~ "Greece",
    pais == "BG" ~ "Bulgaria",
    TRUE ~ pais # Mantiene los nombres que no están en la lista
  ))

Tras el tratamiento de estos datos se puede ver de manera gráfica la comparativa entre los países.

grafico_pib <- ggplot(pib_2022_desc, aes(x = reorder(pais, -`2022`), y = `2022`)) +
  geom_bar(stat = "identity", fill = "gold") +
  labs(x = "País", y = "Valor en 2022", title = "PIB por País en 2022") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1)) # Rotar etiquetas para mejor legibilidad

grafico_pib

Una vez observada la primera variable, generamos los datos del valor medio de la resistencia a antimicrobianos en los países europeos para comprobar si existe una relación o, efectivamente no tiene que ver.

Incidencia_enfermedades <-  read_csv("INPUT/DATA/ECDC_encuesta_AMR_incidencia_enfermedades.csv")

incidencia_2022MF <- Incidencia_enfermedades %>%
  filter(Time == 2022) %>%
  filter(Category == 'Male' | Category == 'Female')

incidencia_2022MF$Value <- as.numeric(incidencia_2022MF$Value)

media_poblacion <- incidencia_2022MF %>%
  select(-Unit, -HealthTopic, -Time, -Distribution) %>%
  mutate(grupo = substr(Population, 1, 3))

lista_codigos_paises <- list("AT", "BE", "BG", "CY", "CZ", "DE", "DK", "EE", "EL", "ES", "FI", "FR", "HR", "HU", 
                             "IE", "IT", "LT", "LU", "LV", "MT", "NL", "PL", "PT", "RO", "SE", "SI", "SK") 


# Calcular la media de 'Value' por 'RegionCode' y 'grupo'
otra <- media_poblacion %>%
  select(-Category, -CategoryIndex, -Population) %>%
  arrange(RegionCode, grupo) %>%    # Ordena los datos
  group_by(RegionCode, grupo) %>%   # Agrupa por RegionCode y grupo
  summarise(mean_value = mean(Value, na.rm = TRUE))%>%
  filter(RegionCode %in% lista_codigos_paises)  # Calcula la media en cada grupo y desagrupa


otra <- media_poblacion %>%
  select(-Category, -CategoryIndex, -Population) %>%
  arrange(RegionCode, grupo) %>%    
  group_by(RegionCode, grupo) %>%   
  summarise(mean_value = mean(Value, na.rm = TRUE))  

otra$mean_value[is.nan(otra$mean_value)] <- 0

media_region <- otra %>%
  group_by(RegionCode) %>%
  summarise(mean_value_region = mean(mean_value, na.rm = TRUE))  


boxplot <- ggplot(otra, aes(x = grupo, y = mean_value)) +
  geom_boxplot() +
  labs(title = "Distribution of Mean Incidence by Bacteria Group",
       x = "Bacteria Group", y = "Mean Incidence") +
  theme_minimal()


incidencia_keyed <- highlight_key(otra, ~RegionCode)

scatter_plot <- ggplot(incidencia_keyed, aes(x = grupo, y = mean_value, color = RegionCode, text = paste("País:", RegionCode))) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE) +
  labs(title = "Population vs Mean Incidence Value by Region",
       x = "Population", y = "Mean Incidence Value") +
  theme_classic()

interactive_scatter_plot <- ggplotly(scatter_plot, tooltip = "text") %>%
  highlight(on = "plotly_click", off = "plotly_doubleclick", color = "red", opacityDim = 0.2)

interactive_scatter_plot

Este código cargado en el entorno, es capaz de generar gráficos que permitan ver la incidencia de las distintas bacterias por cada país y para ver cuál es la bacteria que más afecta en toda la Unión Europea.

boxplot_in <- ggplot(otra, aes(x = grupo, y = mean_value)) +
  geom_boxplot() +
  labs(title = "Distribution of Mean Incidence by Bacteria Group",
       x = "Bacteria Group", y = "Mean Incidence") +
  theme_minimal()
boxplot_in

incidencia_keyed <- highlight_key(otra, ~RegionCode)

scatter_plot <- ggplot(incidencia_keyed, aes(x = grupo, y = mean_value, color = RegionCode, text = paste("País:", RegionCode))) +
  geom_point() +
  geom_smooth(method = "lm", se = FALSE) +
  labs(title = "Population vs Mean Incidence Value by Region",
       x = "Population", y = "Mean Incidence Value") +
  theme_classic()

interactive_scatter_plot <- ggplotly(scatter_plot, tooltip = "text") %>%
  highlight(on = "plotly_click", off = "plotly_doubleclick", color = "red", opacityDim = 0.2)

## `geom_smooth()` using formula = 'y ~ x'

interactive_scatter_plot

Pregunta 3: Consumo y resistencia en sectores específicos en la UE.

Explorar cómo el consumo de antibióticos en sectores como el personal (comunitario), hospitalario y ganadero se relaciona con los niveles de resistencia a antibióticos en diferentes países de la UE.

Uno de los factores a estudiar era el consumo/administración de este tipo de los medicamentos por parte de los hospitales y la ciudadanía de manera individual, y eso se ve reflejado en los siguientes datos extraídos de la encuesta europea anual.

library(tidyverse)
library(rjson)
library(dplyr)
library(tidyr)
library(ggplot2)


DDD_Europa_Json <- fromJSON(file = "INPUT/DATA/DDD_1000_habitantes_paises.JSON")

DDD_Europa_df <- do.call(rbind, lapply(DDD_Europa_Json, function(x) {
  data.frame(Country = x$Country, 
             DDD_per_1000_inhabitants_per_day = as.numeric(x$`DDD per 1000 inhabitants per day`))
}))

DDD_Europa_df <- DDD_Europa_df %>%
  mutate(Country = case_when(
    Country == "Slovakia" ~ "SK",
    Country == "Belgium" ~ "BE",
    Country == "Cyprus" ~ "CY",
    Country == "Greece" ~ "EL",
    Country == "Romania" ~ "RO",
    Country == "Bulgaria" ~ "BG",
    Country == "France" ~ "FR",
    Country == "Malta" ~ "MT",
    Country == "Poland" ~ "PL",
    Country == "Spain" ~ "ES",
    Country == "Ireland" ~ "IE",
    Country == "Italy" ~ "IT",
    Country == "Luxembourg" ~ "LU",
    Country == "Portugal" ~ "PT",
    Country == "Czechia" ~ "CZ",
    Country == "Finland" ~ "FI",
    Country == "Austria" ~ "AT",
    Country == "Norway" ~ "DE",
    Country == "Denmark" ~ "DK",
    Country == "Estonia" ~ "EE",
    Country == "Hungary" ~ "HU",
    Country == "Croatia" ~ "HR",
    Country == "Lithuania" ~ "LT",
    Country == "Latvia" ~ "LV",
    Country == "Netherlands" ~ "NL",
    Country == "Iceland" ~ "SE",
    Country == "Slovenia" ~ "SI",
  ))%>%
  mutate(DDD_per_1000_inhabitants_per_day = DDD_per_1000_inhabitants_per_day/10)

Al tratarse del número de dosis estándar consumidas diariamente por cada 1000 personas, se divide entre 10 para que se convierta en el número de dosis diarias por cada 100 habitantes.

Si se considera que cada dosis se consume por una persona diferente, el resultado de la segunda columna significa que el x% de la población de cada país consume antibiótico de forma diaria en el sector hospitalario y comunitario.

Gráfico de barras de la Dosis de antibióticos media Diaria Definida de los países ordenados de mayor a menor consumo

grafico_DDD <- ggplot(DDD_Europa_df, aes(x = reorder(Country, -DDD_per_1000_inhabitants_per_day), 
                          y = DDD_per_1000_inhabitants_per_day)) +
  geom_bar(stat = "identity", fill = "turquoise") +
  labs(title = "DDD por 1000 habitantes por día en Europa", 
       x = "País", 
       y = "DDD por 1000 habitantes por día") +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))
grafico_DDD

Se puede analizar de forma individual, para determinar qué país consume más antibióticos de forma registrada por el sistema de salud, o por el contrario, se puede comparar con otros factores que se están estudiando como posibles causantes de la resistencia a los antibióticos para encontrar relaciones con el porcentaje de la población de los diferentes países de la UE que han dado positivo en resistencia a los antimicrobianos.

A continuación se pretende estudiar el consumo de antibióticos por parte de los animales de granja suministrados por los profesionales del sector, que en ocasiones resultan en sobremedicación. Este problema se suele dar como resultado de un tratamiento (no necesario) cuyo objetivo es tener animales sanos y sin ninguna enfermedad para poder venderlos y que el producto cárnico no se desperdicie sin perder dinero.

Para abrir los datos relativos a la ganadería, se presenta una función que permite abrir una carpeta con numerosos archivos de extensión .json con la misma estructura y leerlos.

library(httr)
library(tidyverse)
library(jsonlite)
library(readxl)
library(rjson)

obtener_nombre<-function(carpeta){
  archivos <- as.list(list.files(path =carpeta))
  lista_nombres<-list()
  for(i in 1:length(archivos)){
    posicion1<-regexpr("_", archivos[[i]])
    posicion2<-regexpr("\\.", archivos[[i]])
    subcadena<-substr(archivos[[i]], posicion1+1, posicion2-1)
    lista_nombres[[i]]<-subcadena
  }
  
  return(lista_nombres)
  
}


obtener_archivo<-function(direccion){
  lista_paises<-obtener_nombre(direccion)
  lista_enlace<-list()
  direccion_archivos<-list()
  
  for(i in lista_paises){
    cada_pais<-paste0("AMR_",i,".json")
    lista_enlace[i]<-cada_pais
  }
  for(i in lista_enlace){
    
    cada_archivo<-paste0(direccion,"/",i)
    direccion_archivos[i]<-cada_archivo
  }
  
  
  for(i in direccion_archivos){
    pais<-fromJSON(file= i)
    enlace<-pais$links$archive
    respuesta_archivo <- GET(enlace)# Hacer la solicitud HTTP para descargar el archivo
    nombre_archivo<-basename(enlace)#Extrae el nombre del archivo de la URL
    
    if (status_code(respuesta_archivo) == 200) {# Verificar si la solicitud fue exitosa (código 200, código estándar HTTP que significa "OK")
      # Guardar el archivo ZIP localmente en formato binario
      writeBin(content(respuesta_archivo, "raw"), nombre_archivo)
      print("Archivo ZIP descargado correctamente.")
      unzip(nombre_archivo, exdir = "carpeta_destino", overwrite = TRUE)
    } else {
      print(paste("Error al descargar el archivo. Código de respuesta:", status_code(respuesta_archivo)))
    }
    
   
  }
  
}
obtener_archivo("INPUT/DATA/Resistecia_Antibioticos_UE")



leer_archivo <- function(carpeta) {
  carpeta_destino <- carpeta
  archivos_zip <- list.files(carpeta_destino, pattern = "\\.zip$", full.names = TRUE)
  
  # Iterar sobre los archivos .zip y descomprimirlos
  for (archivo in archivos_zip) {
    # Descomprimir el archivo .zip
    archivos_extraidos <- unzip(archivo, exdir = carpeta_destino, overwrite = TRUE)
    print(paste("Descomprimido:", archivo))  # Imprimir cada archivo que se descomprime
    
    # Filtrar el archivo .xlsx entre los extraídos
    archivo_xlsx <- archivos_extraidos[grepl("\\.xlsx$", archivos_extraidos)]#Aquí, grepl() busca archivos cuyos nombres terminen con .xlsx 
    
    # Verificar si hay algún archivo .xlsx descomprimido
    if (length(archivo_xlsx) > 0) {
      # Leer el archivo Excel como dataframe
      datos_xlsx <- read_excel(archivo_xlsx[1])  # Leer el primer archivo .xlsx encontrado
      
      # Asignar el dataframe al Global Environment usando el nombre del archivo como variable
      nombre_variable <- make.names(basename(archivo_xlsx[1]))  # Crear un nombre de variable válido
      assign(nombre_variable, datos_xlsx, envir = .GlobalEnv)  # Asignar el dataframe al Global Environment
      
    } 
  }
}


leer_archivo("carpeta_destino")

Juntar los dataframes que tenemos en el enviroment 1. Listar los nombres de todos los dataframes que terminan en “.xlsx” (ajusta si es necesario) 2. Convertir los nombres a una lista de dataframes usando mget() 3. Unir todos los dataframes en uno solo usando bind_rows

nombres_dataframes <- ls(pattern = "_AMR_PUB\\.xlsx$")


lista_dataframes <- mget(nombres_dataframes)

 
df_combinado <- bind_rows(lista_dataframes, .id = "origen")

Carga del global enviroment El entorno en el que guardan todos los dataframes de cada pais resultantes de la función leer_archivo, y el dataframe df_combinado del codigo anterior se ha cargado anteriormente ya que incluye otros datos que ya se han empleado.

Modificación y obtencion del dataframe con los datos seleccionamos las columnas que vamos a necesitar y filtramos en la columna de las bacterias solo las patogénicas.

Cambiamos los nombres de las columnas y lo asignamos al dataframe

paises_UE_df<-df_combinado%>%
  select(rep_Country_name,rep_Country_code,zoonosis_name,matrix_name,totUnitsTested,totUnitsPositive,sampUnitType_name,sampType_name,MIC_name,CUTOFFVALUE)%>%
  mutate(zoonosis_name = sub(".*", "", zoonosis_name))%>% # Extraer solo la primera palabra
  filter(zoonosis_name != "Escherichia coli, non-pathogenic, unspecified")


nuevos_nombres <- c("NombrePais", "Codigo", "zoonosis_name","OrigenMuestra", "TotalMuestras","MuestraPositiva","Tipo_Unidad_Muestra","TipoMuestra","MIC_name","ValorCorte")  # Modifica según el número de columnas

colnames(paises_UE_df) <- nuevos_nombres

Mapa Interactivo de positivos en ganaderia por paises de la UE

library(sf)
library(ggplot2)
library(plotly)
library(dplyr)
library(viridis)
paises_UE <- c(
  "Cyprus", "France", "Lithuania", "Czechia", "Germany", 
  "Estonia", "Latvia", "Sweden", "Finland", "Luxembourg", 
  "Belgium", "Spain", "Denmark", "Romania", "Hungary", 
  "Slovakia", "Poland", "Ireland", "Greece", "Austria", 
  "Italy", "Netherlands", "Croatia", "Slovenia", "Bulgaria", 
  "Portugal", "Malta"
)



mapa_mudo <- st_read("INPUT/DATA/mapaMundi")  

## Reading layer `ne_10m_admin_0_countries' from data source 
##   `C:\Users\usuario\seminario_fuentes\INPUT\DATA\mapaMundi' using driver `ESRI Shapefile'
## Simple feature collection with 258 features and 168 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -180 ymin: -90 xmax: 180 ymax: 83.6341
## Geodetic CRS:  WGS 84

mapa_mundo_europa <- mapa_mudo %>% 
  dplyr::filter(NAME %in% paises_UE)

mapa_mundo_europa <- mapa_mudo %>% dplyr::filter(NAME %in% paises_UE)


positivos_por_ciudad <- paises_UE_df %>%
  dplyr::filter(NombrePais != "United Kingdom (Northern Ireland)")%>%
  group_by(NombrePais) %>%
  dplyr::summarize(total_pruebas = sum(TotalMuestras, na.rm = TRUE),
            total_positivos= sum(MuestraPositiva, na.rm = TRUE)) %>%
  mutate(ratio_positivo = (total_positivos / total_pruebas) * 100)

# Unir datos de positividad al mapa
mapa_mundo_europa$NAME <- as.character(mapa_mundo_europa$NAME)
positivos_por_ciudad$NombrePais <- as.character(positivos_por_ciudad$NombrePais)

# Realizar el join usando las columnas correctas
mapa_mundo_europa <- mapa_mundo_europa %>% 
  left_join(positivos_por_ciudad, by = c("NAME" = "NombrePais"))



mapa <- ggplot(mapa_mundo_europa) +
  geom_sf(aes(fill = ratio_positivo)) +
  scale_fill_viridis_c(option = "plasma", na.value = "gray") +
  labs(title = "Tasa de Positividad por País en Europa",
       fill = "Tasa de Positividad (%)") +
  coord_sf(xlim = c(-30, 50), ylim = c(35, 72), expand = FALSE) +  # Ajustar límites para hacer zoom en Europa
  theme_minimal() +
  theme(plot.title = element_text(hjust = 0.5))


# Convertir el gráfico a interactivo con plotly
mapa_interactivo <- ggplotly(mapa)

mapa_interactivo